Glacial and erosional contributions to Late Quaternary uplift of the European Alps (GEOLQUEA)

Isostatic adjustments of the Earth’s surface to changes in water, ice, and sediment loading are important contributions to present-day uplift/subsidence rates in many regions on Earth. In the absence of significant horizontal tectonic shortening in the central and western parts of the European Alps, uplift rates larger than 2 mm/yr are difficult to explain by geodynamic processes and have been a matter of debate for many decades. Here we examine the likely contribution of glacial isostatic adjustment in the European Alps in response to changes in ice loading using state of the art ice flow and lithospheric numerical modeling. In contrast to a similar previous approach (Mey et al., 2016), we employ a transient ice sheet model over the last glacial cycle (100 kyr) in combination with a spherical viscoelastic solid earth model. We present ice model results using the Instructed Glacier Model (Jouvet et al., 2021), in which we tested the effect of spatial resolution on the growth and extent of the European ice cap. We found significant differences using a model resolution of 200 m compared to a resolution of 2000 m, which is commonly used in large-scale glacier modeling studies. These differences result in near-steady state volumetric differences at the maximum ice extent of +13% for the high compared to the low-resolution model. In addition, we observed periods of marked ice growth that initiated at significantly different times for the different resolution models. Therefore, we conclude that a realistic ice loading history requires a sufficiently high spatial resolution, which is significantly higher than used in previous models. Based on the modeled ice loading histories, we used the lithosphere and mantle model VILMA (Klemann et al., 2008, J. Geodyn.) to predict the vertical land motion. These estimates are based on a global 60 km thick elastic lithosphere, followed by a 200 km thick viscous layer with a viscosity of 1020 Pa s, which increases to 5 x 1020 Pa s down to 670 km depth, and 3.16 x 1021 Pa s to the core mantle boundary. Preliminary results indicate similar first-order lithospheric responses, with spatiotemporal differences in the magnitude of postglacial response. We hope to present more results based on further ice models that are forced by a more realistic climate history

A visual analytics tool to validate simulation models against collected data. V. 1.0.0

The validation of a simulation model is a crucial task in model development. It involves the comparison of simulation data to observation data and the identification of suitable model parameters. SLIVISU is a Visual Analytics framework that enables geoscientists to perform these tasks for observation data that is sparse and uncertain. Primarily, SLIVISU was designed to evaluate sea level indicators, which are geological or archaeological samples supporting the reconstruction of former sea level over the last ten thousands of years and are compiled in a postgreSQL database system. At the same time, the software aims at supporting the validation of numerical sea-level reconstructions against this data by means of visual analytics

Dependence of late glacial sea-level predictions on 3D Earth structure

Glacial isostatic adjustment is dominated by Earth rheology resulting in a variability of relative sea-level (RSL) predictions of more than 100 meters during the last glacial cycle. Seismic tomography models reveal significant lateral variations in seismic wavespeed, most likely corresponding to variations in temperature and hence viscosity. Therefore, the replacement of 1D Earth structures by a 3D Earth structure is an essential part of recent research to reveal the impact of lateral viscosity contrasts and to achieve a more consistent view on solid-Earth dynamics. Here, we apply the VIscoelastic Lithosphere and MAntle model VILMA to predict RSL during the last deglaciation. We create an ensemble of geodynamically constrained 3D Earth structures which is based on seismic tomography models while considering a range of conversion factors to transfer seismic velocity variations into viscosity variations. For a number of globally distributed sites, we discuss the resulting variability in RSL predictions, compare this with regionally optimized 1D Earth structures, and validate the model results with relative sea-level data (sea-level indicators). This study is part of the German Climate Modeling initiative PalMod aiming the modeling of the last glacial cycle under consideration of a coupled Earth system model, i.e. including feedbacks between ice-sheets and the solid Earth

Data assimilation for a visco-elastic Earth deformation model

We present a data assimilation algorithm for the time-domain spectral-finite element code VILMA. We consider a 1D earth structure and a prescribed glaciation history ICE5G for the external mass load forcing. We use the Parallel Data Assimilation Framework (PDAF) to assimilate sea level data into the model in order to obtain better estimates of the viscosity structure of mantle and lithosphere. For this purpose, we apply a particle filter in which an ensemble of models is propagated in time, starting shortly before the last glacial maximum. At epochs when observations are available, each particle's performance is estimated and they are resampled based on their performance to form a new ensemble that better resembles the true viscosity distribution. In a proof of concept we show that with this method it is possible to reconstruct a synthetic viscosity distribution from which synthetic data were constructed. In a second step, paleo sea level data are used to infer an optimised 1D viscosity distribution

An approach for constraining mantle viscosities through assimilation of palaeo sea level data into a glacial isostatic adjustment model

Glacial isostatic adjustment is largely governed by the rheological properties of the Earth's mantle. Large mass redistributions in the ocean–cryosphere system and the subsequent response of the viscoelastic Earth have led to dramatic sea level changes in the past. This process is ongoing, and in order to understand and predict current and future sea level changes, the knowledge of mantle properties such as viscosity is essential. In this study, we present a method to obtain estimates of mantle viscosities by the assimilation of relative sea level rates of change into a viscoelastic model of the lithosphere and mantle. We set up a particle filter with probabilistic resampling. In an identical twin experiment, we show that mantle viscosities can be recovered in a glacial isostatic adjustment model of a simple three-layer Earth structure consisting of an elastic lithosphere and two mantle layers of different viscosity. We investigate the ensemble behaviour on different parameters in the following three set-ups: (1) global observations data set since last glacial maximum with different ensemble initialisations and observation uncertainties, (2) regional observations from Fennoscandia or Laurentide/Greenland only, and (3) limiting the observation period to 10 ka until the present. We show that the recovery is successful in all cases if the target parameter values are properly sampled by the initial ensemble probability distribution. This even includes cases in which the target viscosity values are located far in the tail of the initial ensemble probability distribution. Experiments show that the method is successful if enough near-field observations are available. This makes it work best for a period after substantial deglaciation until the present when the number of sea level indicators is relatively high

Validating GIA models based on an ensemble of 3D Earth structures with present-day GPS uplift rates

Glacial-isostatic adjustment (GIA) models simulate the viscoelastic response of the solid earth due to loading. During the last glacial maximum, large areas in the northern and southern hemisphere were covered by km-thick ice sheets. Although most of the ice has been melted already 8,000 year ago, the time-delayed response of the viscoelastic earth is still a significant contribution to present-day uplift rates. The implementation of GIA models in global climate models is an essential part of the current research. Hereby, the choice of an appropriate earth structure in the GIA model plays an important role and has to be constrained by observational data. Here, we apply present-day uplift data to constrain a set of GIA models that differ in 3D earth structure. To this end, these different GIA models are validated against GPS uplift rates provided by Schumacher et al. (2019). The GPS stations are globally distributed and not necessarily clustered in regions with strong GIA signal. For validation, regions with the largest gradient present in the GIA signal are most crucial. Thus, we use a weighting scheme, where those GPS stations get a higher weight that are less correlated to all other stations. Additionally, uncertainties in the GPS rates appear due to the length of the GPS time series and due to station specifics such as the used GPS receiver, and are provided together with the rates as standard deviations. Thence, the weighting used for the validation is the sum of the correlation derived weights and the uncertainty derived weights. With this weighting in place, different GIA models can be validated against present day uplift rates by means of root mean square errors or mean absolute error

The Influence of Sediments, Lithosphere and Upper Mantle (Anelastic) With Lateral Heterogeneity on Ocean Tide Loading and Ocean Tide Dynamics

Ocean tide loading (OTL) and ocean tide dynamics (OTD) are known to be affected by Earth's internal structures, with the latter being affected by the self-attraction and loading (SAL) potential. Combining the 3D earth models Lyon and LITHO1.0, we construct a hybrid model to quantify the coupled effect of sediments, oceanic and continental lithosphere, and anelastic upper mantle on OTL and OTD. Compared to PREM, this more realistic 3D model produces significantly larger vertical OTL displacement by up to 3.9, 2.6, and 0.1 mm for the M2, K1, and Mf OTL, respectively. Moreover, it shows a smaller vector difference of 0.1 mm and a smaller amplitude difference of 0.2 mm than PREM with OTL observations at 663 Global Navigation Satellite System stations, a confirmation of the cumulative effect due to these earth features. On the other hand, we find a resonant impact of wider extent and larger magnitude on OTD, especially for the M2 and K1 tides. Specifically, this impact is concentrated in the ranges 0–6 mm and 0–1.5 mm for M2 and K1, respectively, which is considerably larger than the impact on SAL (mostly in the ranges 0–2 mm and 0–1.0 mm, respectively). Since the effect on vertical displacement is at a similar level compared to the accuracy of modern data-constrained ocean tide models that require correction of the geocentric tide by loading induced vertical displacements, we regard its consideration to be potentially beneficial in OTD modeling

Gravitationally Consistent Mean Barystatic Sea Level Rise From Leakage‐Corrected Monthly GRACE Data

Gravitationally consistent solutions of the Sea Level Equation from leakage‐corrected monthly‐mean GFZ RL06 Gravity Recovery and Climate Experiment (GRACE) and GRACE Follow‐On (GRACE‐FO) Stokes coefficients reveal that barystatic sea level averaged over the whole global ocean was rising by 1.72 mm a−1 during the period April 2002 until August 2016. This rate refers to a truely global ocean averaging domain that includes all polar and semienclosed seas. The result corresponds to 2.02 mm a−1 mean barystatic sea level rise in the open ocean with a 1,000 km coastal buffer zone as obtained from a direct spatial integration of monthly GRACE data. The bias of +0.3 mm a−1 is caused by below‐average barystatic sea level rise in close proximity to coastal mass losses induced by the smaller gravitational attraction of the remaining continental ice and water masses. Alternative spherical harmonics solutions from CSR, JPL, and TU Graz reveal open‐ocean rates between 1.94 and 2.08 mm a−1, thereby demonstrating that systematic differences among the processing centers are much reduced in the latest release. We introduce in this paper a new method to approximate spatial leakage from the differences of two differently filtered global gravity fields. A globally constant and time‐invariant scale factor required to obtain full leakage from those filter differences is found to be 3.9 for GFZ RL06 when filtered with DDK3, and lies between 3.9 and 4.4 for other processing centers. Spatial leakage is estimated for every month in terms of global grids, thereby providing also valuable information of intrabasin leakage that is potentially relevant for hydrologic and hydrometeorological applications

Joint inversion estimate of regional glacial isostatic adjustment in Antarctica considering a lateral varying Earth structure (ESA STSE Project REGINA)

A major uncertainty in determining the mass balance of the Antarctic ice sheet from measurements of satellite gravimetry, and to a lesser extent satellite altimetry, is the poorly known correction for the ongoing deformation of the solid Earth caused by glacial isostatic adjustment (GIA). Although much progress has been made in consistently modelling the ice-sheet evolution throughout the last glacial cycle, as well as the induced bedrock deformation caused by these load changes, forward models of GIA remain ambiguous due to the lack of observational constraints on the ice sheet's past extent and thickness and mantle rheology beneath the continent. As an alternative to forward modelling GIA, we estimate GIA from multiple space-geodetic observations: GRACE, Envisat/ICESat and GPS. Making use of the different sensitivities of the respective satellite observations to current and past surface mass (ice mass) change and solid Earth processes, we estimate GIA based on viscoelastic response functions to disc load forcing. We calculate and distribute the viscoelastic response functions according to estimates of the variability of lithosphere thickness and mantle viscosity in Antarctica. We compare our GIA estimate with published GIA corrections and evaluate its impact in determining the ice mass balance in Antarctica from GRACE and satellite altimetry. Particular focus is applied to the Amundsen Sea Sector in West Antarctica, where uplift rates of several cm/yr have been measured by GPS. We show that most of this uplift is caused by the rapid viscoelastic response to recent ice-load changes, enabled by the presence of a low-viscosity upper mantle in West Antarctica. This paper presents the second and final contribution summarizing the work carried out within a European Space Agency funded study, REGINA, (www.regina-science.eu)

An approach for constraining mantle viscosities through assimilation of palaeo sea level data into a glacial isostatic adjustment model

Glacial isostatic adjustment is largely governed by the rheological properties of the Earth's mantle. Large mass redistributions in the ocean–cryosphere system and the subsequent response of the viscoelastic Earth have led to dramatic sea level changes in the past. This process is ongoing, and in order to understand and predict current and future sea level changes, the knowledge of mantle properties such as viscosity is essential. In this study, we present a method to obtain estimates of mantle viscosities by the assimilation of relative sea level rates of change into a viscoelastic model of the lithosphere and mantle. We set up a particle filter with probabilistic resampling. In an identical twin experiment, we show that mantle viscosities can be recovered in a glacial isostatic adjustment model of a simple three-layer Earth structure consisting of an elastic lithosphere and two mantle layers of different viscosity. We investigate the ensemble behaviour on different parameters in the following three set-ups: (1) global observations data set since last glacial maximum with different ensemble initialisations and observation uncertainties, (2) regional observations from Fennoscandia or Laurentide/Greenland only, and (3) limiting the observation period to 10 ka until the present. We show that the recovery is successful in all cases if the target parameter values are properly sampled by the initial ensemble probability distribution. This even includes cases in which the target viscosity values are located far in the tail of the initial ensemble probability distribution. Experiments show that the method is successful if enough near-field observations are available. This makes it work best for a period after substantial deglaciation until the present when the number of sea level indicators is relatively high